COMMENT ⊗   VALID 00006 PAGES
C REC  PAGE   DESCRIPTION
C00001 00001
C00002 00002	_Dynamic Programming - Exercise (Recommended)_
C00006 00003	_Recursion_
C00017 00004		ERASEAT(3,7)
C00020 00005	Rules of good program practice
C00022 00006	_Tournament Sorting_
C00031 ENDMK
C⊗;
_Dynamic Programming - Exercise (Recommended)_

The latest electronic gambling game in Las Vegas is a coin flipping machine.
Put your silver dollar in the slot; the machine flips its coin and pays you
two dollars if it comes up heads.  Because of very undisciplined programming
at M-A Mechanics, Inc., makers of the machine, each machine has a definite
probability of tossing heads, but the probability varies widely from machine
to machine; in fact, a survey of the machines suggests that the probabilities
are uniformly distributed between 0 and 1.  A statistician remarks that if
you have played several times against an otherwise unknown machine of this 
type, the above implies that after H heads and T tails, the chance of heads
next time is (H+1)/(H+T+2).

I am changing planes in Las Vegas, and I have time to bet 100 coins on the
one coin-flipping machine in the airport.  It seems intuitively right to keep
betting as long as I'm winning, and to quit if my losses get too big.  But
how big is too big?  If the machine has tossed H heads and T tails, should I
bet or quit?  Print a table that gives me the right strategy.  Don't assume 
that anything about this problem is obvious.  After winning 4 and losing 5 coins,
I should definitely play one more time, for example.

Please hurry; my plane leaves in twenty minutes.

_Exercise (Harder)_

The firm of Muscle Victimless Enterprises has just delivered a second coin
flipping machine to the airport lounge.  Now what should I do?  I only have
time to play twenty coins now.  To keep the size of the printout down, just
print situations, if there are any, in which I should not play the machine
currently more likely to throw heads.  For example, machine 1 has had heads
5 times, tails 4; machine 2 has had heads and tails once each.  If I try
machine 2 and it throws heads, it will start looking much better than machine 1;
maybe I should give it a chance.

_Exercise - Dynamic Programming_

The thirtieth Fibonnacci number is __________; to compute it iteratively takes
thirty iterations.  This recursive function can compute it as FIB(30):

	FUNCTION FIB(I:INTEGER): INTEGER
		BEGIN
		IF I<2 THEN FIB:=I
		ELSE FIB:=FIB(I-1)+FIB(I-2)
		END;

how many (recursive) calls to the procedure are needed to do the computation?

_Recursion_

How do you take the derivative of a long formula?  I break it into several
shorter formulas, take their derivatives, and combine the results.  If the 
formula is  F=G/H,  I get G' and H' in separate calculations, then combine
them as  F'=(G'H-H'G)/H↑2.   The process, superficially, looks circular;
before I can find one derivative, I have to find another.  In fact the
subordinate derivations get easier and easier, acting on shorter formulae;
eventually all the subproblems reduce to finding derivatives of variables
or constants, which can be done directly.

An algorithm which solves problems by taking time out along the way to solve
simpler problems of the same kind, using the same methods, is called _recursive_.
To execute such an algorithm, you must keep track separately of variables for
each subordinate problem, and of where you are in the instructions for each
subordinate problem.  For example, you could write each subordinate problem
on a separate sheet of paper, and cross reference the sheets to allow returning
from a subordinate problem to its parent problem.  As we shall see, Pascal
procedures and functions can be used to carry out recursive algorithms.

*******Insert adaptive integration algorithm.

Another recursive algorithm can be used to allow input of real numbers as
arbitrary arithmetic formulae built up from real constants.  For example,
if a program requires an input of a temperature measured in Celsius units,
and I want to enter 98.6 \degrees Fahrenheit, it would be convenient to type
((5*(98.6-32)/9) rather than hand-calculating the converted value.  The
algorithm is this:

If the first character is '(', the input must be of the form (F O G), where
F and G are formulae (possibly simply numbers) and O is one of +, -, *, /.
Using the same algorithm, read F, getting a numerical value Y; read the
operator O and save it; read G, getting a numerical value Z; read the ')';
depending on the operator character O, return  Y+Z, Y-Z, Y*Z, or Y/Z.
If the first character is not (, read a number and return it.

In Pascal, this algorithm is expressed as a procedure:

PROCEDURE READFORM(F:TEXT; VAR X: REAL);
VAR Y,Z: REAL; OP: CHAR;
BEGIN
(IF F↑<>'' THEN
	READ(F,X)
ELSE
	BEGIN
	(GET(F); (*DISCARD'('*)
	READFORM(F,Y);
	READ(F,OP);
	READFORM(F,Z);
	IF F↑<>')' THEN WRITE(error message);
	GET(F);
	CASE OP OF	(*FIX UP SYNTAX*)
		'+': X:=Y+Z;
		'-': X:=Y-Z;
		'*': X:=Y*Z;
		'/': X:=Y/Z END(*CASE*)
	END
END

To see what the procedure does in detail, insert at the beginning and end of the
procedure body the commands
	WRITELN('STARTING READFORM') and
	WRITELN('LEAVING READFORM').
After each read, echo the input; after READFORM(F,Y), say, include the command
	WRITELN('RESULT Y OF RECURSION', Y).
I recommend reading the output of executing this expanded procedure on several
formulae.

Next we look at how a recursive procedure can be executed by hand.  Suppose
we have a character array, C[1..N,1..N], which contains a number of non-blank
regions, or blobs.  We say that two non-blank characters must belong to the
same blob if they are adjacent horizontally or vertically.  In this diagram,

		A B
		C  E F
		D  G H

ABCD is one blob, EFGH another.  As part of a procedure to count the blobs in
a region, we want a procedure to check if there is a non-blank character in a
specified position, and if so, to erase the entire blob to which that character
belongs.

The algorithm to erase a blob at array coordinates X,Y is readily expressed
recursively:

	If C[X,Y] is blank, do nothing.
	If C[X,Y] is non-blank, make it blank and apply the same algorithm
	to the four adjacent points.

We assume there is a border of blanks around the picture to avoid going out
of bounds.

To show the algorithm is a correct one, we show it always halts, that it erases
everything connected to the original point, and that it erases nothing else.

The body of the algorithm has no iterations, so the algorithm could go on forever
only by creating infinitely many subproblems.  At each point, only the first
application of the algorithm to the point can create subproblems; later ones
find C[X,Y] blank and immediately are finished.  Therefore at most N↑2 of the
applications can create subproblems, and at most 4n↑2 subproblems can be created.

If the point (U,V) is initially in the same blob as (X,Y), there is a non-blank
path between them.  If the execution of the algorithm does not reach (U,V),
there is a first place on the path it fails to reach.  But when it first
reaches the previous point on the path, that point is non-blank, and it creates
a subproblem at all four surrounding points; this would be a contradiction, so
the algorithm reaches all points of the blob.

The algorithm has no way to get outside the set of points consisting of the blob
itself and the blanks immediately adjacent to it, so it erases the blob, the
whole blob, and nothing but the blob.

In Pascal, the algorithm is expressed as:

PROCEDURE ERASEAT(X,Y: INTEGER);
BEGIN
IF C[X,Y]<>' ' THEN
	BEGIN
	C[X,Y]:=' ';
	ERASEAT(X+1,Y);		(*1*)
	ERASEAT(X-1,Y);		(*2*)
	ERASEAT(X,Y+1);		(*3*)
	ERASEAT(X,Y-1)
	END
END

To execute by hand the application of the procedure by the call

	ERASEAT(3,7) on the array shown below,	

	Y = 1 2 3 4 5 6 7 8 9 10
    X =	1
	2
	3		* *
	4		*
	5		*
	6
	7
	8
	9
       10

we start by creating a contour box

	ERASEAT(3,7)

	X = 3	Y = 7	Return to Main





showing the procedure being executed and the values of its arguments outside;
inside, current values of parameters, and where to return when finished are
recorded.  Executing the procedure body, we set C[3,7] to blank, and make
another call on ERASEAT, starting a new contour box; we now have:


	ERASEAT(3,7)

	X = 3	Y = 7	Return to Main

	ERASEAT(4,7)

	X = 4	Y = 7	Return to 1




There are now two distinct variables, both called X, with different values.
One, with the value 4, is in current use; until we close the inner box, any
reference to X is assumed to mean the inner one, so the next call,
ERASEAT(X+1,Y) means ERASEAT(5,7).  After completion of the second call, shown
by closing its box, reference to X is assumed to mean the outer one, resuming
the previously interrupted procedure call where X is 3.

After several more steps, the contours look like this, leaving out some
redundant information:
	ERASEAT(3,7)

	X 3	Y 7	R Main
	  erases C[3,7]

	X 4	Y 7	R 1
	  erases C[4,7]

	   X 5	Y 7	R 1
	     erases C[5,7]

	     X 6   Y 7   R 1

	     X 4   Y 7   R 2

	     X 5   Y 8   R 3

	     X 5   Y 6   R 4

We return to the level where X = 5, Y = 7, at the point labelled (*4*), and find
it is finished; we close off the bottom of its box,




and return to the level where X = 4, Y = 7, at (*1*)

             X 3   Y 7	 R 2

	     X 4   Y 8	 R 3

	     X 4   Y 6	 R 4




Now we are back at the original call, with X = 3, Y = 7, at (*1*)



	     X 2   Y 7	 R 2

	     X 3   Y 8	 R 3
	      erases C[3,8]

	         X 4   Y 8   R 1

	         X 2   Y 8   R 2

         	 X 3   Y 9   R 3

	         X 3   Y 7   R 4


	     X 3       Y 6   R 4



Now we are back in the main program, having completed the execution of
ERASEAT(3,7).  The variables X and Y no longer exist, as indicated by the
closing of all the boxes containing them.  The blob has been erased.

Naturally, computer execution of recursive procedure is not literally done by
drawing boxes.  It does involve keeping values for several variables of the
same name, and remembering where in the program execution should return after
each procedure call is finished.
Rules of good program practice

Good programmers do some of their calculations by hand.  A computer program is
like a submarine; you know it's there, but you can't see what it's doing.
Until it's too late.

Human beings, unlike computers, have the gift of recognizing patterns in events.
They are also blessed by laziness.  To evaluate the polynomial a+bx+cx↑2+dx↑3+ex↑4
in the obvious way requires seven multiplications and four additions. If you have
to do it by hand a few dozen times, you get strongly motivated to notice a common
factor, and calculate a+x(b+cx+dx↑2+ex↑3) using one fewer multiplication.  From
there, it is not far to Horner's Rule: a+x(b+x(c+x(d+x(e)))), saving four
multiplications.

The next best thing to executing an algorithm by hand is to get as much information
as you can about how a computer executes it.  You can use an array of counters to
find out how often each iteration body and procedure body in your program is
executed.  An old programmer's rule of thumb declares that 80% of the time is
spent in 20% of the program, and 80% of _that_ time is spent in 20% of the 20%.
Find that hyperactive 4% and make it faster; you will speed your entire program
up by almost as much.
_Tournament Sorting_

Some of our example algorithms have shown that searching a table is much faster
if that table is sorted.  How much would you pay for a dictionary or a phone
book if they were not alphabetized?  Here is an algorithm for sorting that is
simple, yet efficient.  For  n  data, make a tree like the one below, with
2n-1 circles:

					1



			2				3



		4		5		6		7



	   8	    9 	  10       11	   12       13 


In the example, n=7.  Put your data in the last  n  circles, which have no other
circles hanging from them,



					1



			2				3



		4		5		6		7
								Sunday


	   8        9	   10       11	   12      13
	   Mon      Tue	   Wed      Thu	   Fri     Sat


Treat the diagram like the pairings for a tournment.  Player 12, Fri, plays
Player 13, Sat.  The winner is Fri, because her name is alphabetically earlier;
Fri's name and number go in circle number 6, and we continue counting down
until we have a winner in circle 1:


					12
				1	Fri



			8				12
		2	Mon			3	Fri



		8		11		12	 	7
	4	Mon	5	Thu	6	Fri	7	Sun



	   8	    9	   10       11 	     12        13
      8	   Mon  9   Tue 10 Wed  11  Thu  12  Fri   13  Sat            


The winner is Fri, player 12.  We write or record the winner, then disqualify Fri
to see who comes in second.  To disqualify Fri, we change circle 12 to ZZZ and
work upward in the diagram, through circles 6,3, and 1, finding new winners:


					8
				1	Mon


							13
			2			3	Sat


					     13
		4		5	6    Sat	7



					     12          13
	8		9    10    11	 12  ZZZ    13   Sat
					

				(others as above)

We keep recording and disqualifying winners until we have counted  n  winners
(or until the winner is ZZZ).

To program this algorithm, we must be able to follow the connections up and down.
The circles are numbered from 1 to 2n-1.  If a circle has number  i < n, the two
below it are numbered  2i  and  2i+1.  If  i > 1, the one above it is numbered
i DIV 2.  In Pascal, the algorithm follows:

	PROGRAM TOURNAMENT (		);
	declarations
	BEGIN
	FOR I:= 1 TO 2*N+1 DO NUMBR[I]:=I;
	FOR I:= N TO 2*N+1 DO READ(PLAYER[I]);
	FOR I:= N-1 DOWN TO 1 DO
		IF PLAYER [2-I]<PLAYER[2*I] THEN
			BEGIN
			NUMBER[I]:=NUMBER[2*I]
			PLAYER[I]:=PLAYER[2*I]
			END
		ELSE
			BEGIN
			-----------------[2*I+1]   see above
			-----------------[2*I+1]
			END;
	WHILE PLAYER[1]<'ZZZ' DO
		BEGIN
		WRITE(PLAYER[1]);
		J:=NUMR[1];
		PLAYER[J]:= 'ZZZ';
		WHILE J>1 DO
			BEGIN
			I:= J DIV 2;
		IF PLAYER[2*I]----------THEN



		ELSE



		END;
		J:= I
		END(*WHILE J*)
	END(*WHILE PLAYER*)
END(*PROGRAM*)

This algorithm, applied to 1000 data, finds the first with 999 comparisons on the
input data.  Thereafter, every other item on the sorted list is produced with at
most ten comparisons, for a total of about 11000.  Most comparably simple sorting
algorithms would require 5000000 comparisons.

The arrays used by the tournament sorting algorithm have an implicit structure of
connections among th array elements, allowing us to get from an array element to
the elements above or below it in the tree.  Such a structure of connections
among data is called, unsurprisingly, a _data structure_.  The most common data
structure of an array in which element  i  is connected to element  i+1.  There
are many other structures, simple and complicated, on which useful algorithms
have been based.

_Successive Approximation_

The square root algorithm described above is an example of a method of computing
which begins with an approximation to the answer, whether good or poor, and 
applies to the approximate answer an operation which is certain to improve it
if it can be improved.  For a problem to which the solution is a real number,
the pattern required is:

   (1) If x is not the solution, f(x) is closer to the solution than  x  is.

   (2) If x is the solution, f(x)=x.

   (3) The function  f  is continuous.

One can show that repeated application of  f  to any starting value gives results
which approach the solution as closely as desired.

If the solution to a problem is an integer, condition (3) is not required.  In
practice, f(x) should be closer to the solution than  x  is by at least a
substantial constant factor; if  f(x)  is less than half as far as  x  from the
solution, thirty or so iterations gives the solution to nine decimal places.
The square root algorithm, with a starting value in error by a factor of 2,
gives the solution, to nine places, in about four iterations.

_Exercise_

If  f  has property 2, and if  |df/dx|<1, show it has properties (1) and (3).


_Exercise_

Find the positive solution to the equation  ax↑3+bx↑2+cx+d=0  by using
f(x)= -(ax↑3+bx↑2+d)/c.  Show that  f  has the desired properties in an
interval that includes the solution.


_Exercise_

Let x↓0 and x↓1 be any two positive numbers, and x↓{i+2}=x↓{i+1}(i≥2).
Show the ratio  x↓{i+1}/x↓i  approaches, as a limit, the positive solution of
the equation  x↑2=x+1.
(Note:  x↓{i+2}/x↓{i+1}=1+1/(x↓{i+1}/x↓i), or  r↓{i+1}=f(r↓i)=1+1/r↓i.


_Exercise_

The function  g(x) has
a root at r.  A common
way to find such a root,
given an estimate x, is
to evaluate the function
and its derivative at x,
and use as the next
estimate the zero of a
linear function tangent to
g at x; this gives
f(x)= x-g(x)/g'(x).
Find some conditions under which this method can be guaranteed to converge.

f' = 1 - (g')↑2 - gg''  =    gg''   ;
	   (g')↑2	   (g')↑2


f' < 1 if gg''<(g')↑2, or  g<(g')↑2/g''.